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Abstract. In this paper we develop the theory of high-order radiation boundary 
conditions for wave propagation problems. In particular, we study the convergence of 
sequences of time-local approximate conditions to the exact boundary condition, and 
subsequently estimate the error in the solutions obtained using these approximations. 
We show that for finite times the Pade approximants proposed by Engquist and Majda 
lead to exponential convergence if the solution is smooth, but that good long-time error 
estimates cannot hold for spatially local conditions. Applications in fluid dynamics are 
also discussed. 
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1. Introduction. Problems in wave propagation are generally posed 
on unbounded domains. Their numerical solution thus requires the intro- 
duction of an artificial boundary and the imposition of radiation boundary 
conditions there. Scores of authors have considered this problem, and a 
number of reasonably accurate procedures have been discovered. Nonethe- 
less, in order to obtain some specified accuracy, it is still generally the 
practice to enlarge the domain - a process which may be inefficient and 
difficult to automate. 

In this work we pursue a different approach - namely to fix the artifi- 
cial boundary and to improve the accuracy by increasing the order of the 
approximate radiation conditions. From a practical point of view, we see 
that these high-order conditions can be easily implemented via the intro- 
duction of auxiliary functions on the boundary. From a theoretical point of 
view, estimates of convergence for fixed boundaries and increasing order are 
needed. We develop such estimates for the wave equation in a half-space 
by first finding a convenient representation of the exact radiation condi- 
tion, which turns out to involve convolution in time with a Bessel kernel. 
Approximate conditions are similarly represented in terms of convolutions, 
and the error then depends on the difference between the exact and ap- 
proximate kernels. Using approximations to an integral representation of 
the exact kernel, convergent time-local approximate conditions are derived. 
These include the spatially local Pade conditions proposed by Engquist and 
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Majda [7,8]. For long time computations, on the other hand, it is shown 
that spatially nonlocal conditions are generally needed. 

Generalizations to other problems are presented, including the lin- 
earized Euler equations as well as the wave equation with circular and 
spherical boundaries. Throughout we indicate some interesting theoretical 
and practical issues which remain unresolved. 

2. The Wave Equation in a Half-Space. 

2.1. Exact Boundary Conditions. We consider: 

(2.1) utt = c 2 V 2 u + /, t > 0, x = (zi,y) € (0,oo) x R n ~ l , 

(2.2) u(x, 0) = g(x), B c u( 0, y, t) = g 0 (y } t). 

We suppose, for some X, 6 > 0, that / = g = 0 and c = 1 for X\ > L — 6. 
Let u(aei, Jb, $) be the Fourier- Laplace transform of u with respect to y and 
i. Then it is easily shown (e.g. [11]) that & satisfies the exact boundary 
condition at X\ = L: 

(2.3) ^ + ( S 1 +|*| 1 ) 1/1 « = 0. 

(The branch of (s 2 + IJb] 2 ) 1 / 2 is chosen so that it is analytic in the right 
half s-plane and has positive real part.) This exact condition is expressed 
in terms of u in the following way: Let T denote Fourier transformation 
with respect to y and JF -1 be its inverse. Let: 

(2.4) £(t) = = — f \/T— "^cosurfdu/. 

As shown in the appendix, 

(2.5) £(a) = (a 2 + l) 1/2 - a. 

Using standard formulas from Laplace transform theory (e.g. [5]) we finally 
have the exact condition at X\ = L: 

(2.6) ^ ^ + F- 1 {\k\ 2 IC(\k\t) * (Tu(x u .*))) = 0. 

(Here, * denotes convolution.) 

2.2. Approximate Conditions. Although it may be possible to di- 
rectly implement (2.6) using FFT’s and fast convolutions, most work has 
been focussed on the development of approximate conditions involving dif- 
ferential operators in time and, usually, space. Local approximations in 
time correspond to rational approximations in s to /C: 

(2.7) (a 2 + |Jfe| 2 ) 1/J - * = |fc| (C* 2 + 1) 1/J - z) « |*|*( 2 ), * = s/|*|. 



3 


We take R to be a rational function of degree (p,p+ 1), that is, 

(2.8) A(z) = deg(P) = p, deg(Q) = p + 1. 

The approximate condition may be directly localized in time by apply- 
ing the operator T~ x Q{\k\~ 1 d / dt)T . However, this leads to differential 
operators of high order as p is increased. To develop a more convenient 
framework for implementation, we make the additional assumption that 
the roots of Q are distinct. Then, R has a partial fraction expansion: 


(2.9) 


Let 

( 2 . 10 ) 


R(z) 


E 


i = i 


a j 

Z-Pj' 


= -SiiiU, 

Z-Pj 


and let |fcj,t) be the Fourier transform of v with respect to y . Here, 
v denotes the approximation to u computed on the bounded domain. We 
finally have the approximate boundary condition: 


( 2 . 11 ) 



= 0 , 


(2.12) (J^- l*|Pj) hj = aj [fc| 2 5. 

The advantage of this formulation is clear: the order is increased simply 
by increasing the number of terms in the sum. From the point of view of 
code development, this is very convenient. 

We note that the conditions above are still nonlocal in space. For 
periodic problems this is no obstacle, as FFT’s can be used. However, the 
nonlocality does preclude their use at more general boundaries. A glance 
at (2.12) reveals the condition for spatial locality: the poles, py, of R{z) 
must come in conjugate, imaginary pairs, or be 0 and R itself must be an 
odd function of z. We may then assume that R has an expansion of the 
form: 


< 2 - i3 > = i 

This leads to the local implementation: 


dv 

dxi 



= 0, 


(2.14) 
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( 2 - 15 ) (|j -^ V Un) 

We would like to impose such a locality condition on <R(s), but it will be 
shown later that such approximations cannot lead to good error estimates 
uniformly in time or in tangential wave number. 

In what follows, we will view 

(2.16) Q{t) = C~ l R{s), 

as an approximation to K(t ). For reference we note that (2.9) corresponds 
to: 

p+i 

(2.17) G(t) = 

i - 1 


while (2.13) implies: 

* 

(2.18) G{t) = ^7; cos fyt. 

J=i 

2.3. Error Estimates. We have seen that it is relatively straightfor- 
ward to implement conditions of increasing order, at least in the half-space 
(or periodic) case. This leads to the question of convergence. Naturally, 
error estimates for approximate boundary conditions have been considered 
(e.g. [2], [8], [13], [20]). However, none of these consider convergence for a 
fixed problem in a fixed domain as the order of the conditions is increased. 
Let e = u — v be the error. Then e satisfies: 

(2.19) e tl = c 2 V 2 e, t > 0, s = (®i,y)G(0,I)x JT" 1 = 

(2.20) «(*, 0) = 0, B o e(0, y , t) = 0, 

(2.21) ^ + |fc| 2 £(l*IO * « = (1*1*) * *• 

The error kernel is given by: 

( 2 . 22 ) £{r)=G{r ) - 

Estimates of e naturally require both the stability and consistency of 
the approximate boundary conditions. Stability is a consequence of the 
uniform Lopatin ski condition: 


(2.23) 


s + {s 2 + \k\ 2 ) l ' 2 + \k\R(z)tQ, 
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for 

(2.24) &($) > 0, k 6 IT- 1 , (5, k ) ^ (0, 0). 

Then we have (e.g. Sakamoto [22, Ch. 3]): 

T T 

(2-25) Jq ||e(-,f)||^r J [n]dt < C 3 ^ \\Eu(L, >t)llfr 0 [ji — >]^i 
where, 

(2.26) fJiu = ^- 1 |Jb| a £(|fc|t)*(^). 


The error may now be bounded in terms of the error in the approximation, 
Q y to /C. In particular, suppose, for some /j > 0 and T > 1: 


(2-27) II^IU.KO.T)] < «T»*. 

Then, by Parseval’s identity and standard estimates for convolutions, 

rp 

J ||Su(X,,t)|||r o[J? .- l] dt = j j^Jk\\£{\k\t)*u{L,k,t)){t)dkdt 

(2-28) < f f l^ ! l 2 | s (^i^>*)| 2 ||^lli,[(o,|t|T)] 

J 0 JR* -1 

< [ T ML, -, 

JO 

Substituting this into (2.25) we finally obtain: 

( 2 * 29 ) ^ € 2 T 2fi C 2 ||u(£> 

This error estimate is best, both from the point of view of long time 
behavior and from the point of view of smoothness required of u, if p = 0. 
We note that such an estimate requires bounds on H^llLiKo.oo)]* This can- 
not be attained for local conditions, as we have seen that they involve 
convolution kernels which are combinations of cos /3 ; *i (2.18), and, hence, 
are not elements of Zi[(0, oo)]. Time uniform estimates could be obtained 
using spatially nonlocal conditions, however. Some discussion of long-time 
behavior of spatially nonlocal boundary conditions appears in [6,17]. In 
[15] we construct conditions using Laguerre and exponential expansions. 
Although the conditions so derived do lead to estimates with fx = 0, con- 
vergence as the order of approximation was increased was slow at best. 
Below we introduce a new nonlocal approximation based on the direct ap- 
proximation to an integral representation for /C(t). 

Our point of view leads to an interesting, and to our knowledge un- 
solved, problem in approximation theory. Define 7i p to be the set of all 
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real functions in £i[(0, oo)] whose Laplace transform is a rational function 
of degree (p,p + 1). More directly, a function in Tip takes the form: 

(2.30) G(t) = cos (fit + d>j) + £ 

i= i J =1 

where 7 and Kj are positive and Mj(t) and Nj(t) are polynomials. Here 
p is given by: 


Pi 

(2.31) P + 1 = 2 53(deg(Mj ) + 1) 


P 3 

1 _L \ ( rlo tr( AT_- ^ 


J=1 


J = 1 


Problem A: For fixed p characterize and find an algorithm to compute the 
best Li[(0, oo)] approximation, Q £ Tip, to 1C or to other kernels. Estimate 
the behavior of the error as p is increased. 

The solution of this problem would provide us with optimal approxi- 
mations to convolutions via the solution of differential equations. We note 
that for bounded intervals and sums excluding trigonometric terms (i.e. 
Pi = 0), a theory does exist. (See Braess [3, Ch. 6].) 

2.4. Methods Derived via Quadrature. In general, approximate 
boundary conditions have been derived either by direct approximation to 
the symbol (e.g. [24]) or through the use of far-field asymptotics (e.g. 

[2,16]). Here we show how a class of convergent local (in space and time) 
approximate conditions may be derived by approximating the integral rep- 
resentation of the exact kernel: 


(2.32) 


= — f yf 1 — w 2 cos urtdw. 

t W-1 


The simplest example is the trapezoid rule: 


(2.33) 


- 

t 



1 — Wj COS Wjt, 


Wj = -1 + 


JL 

9+1* 


After Laplace transformation we find that: 


(2.34) 


R(z) 



\f l - w i z 

2 2 + W] ‘ 


(Of course, in implementations of the condition the number of terms in the 
sum can be halved using the evenness in Wj of the integrand.) 

The well-posedness of these approximations follows directly from a 
check on the Lopatinski condition. (It also follows from general results 
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for local conditions in [12,23].) In order to derive error estimates, however, 
we must explicitly bound the stability constant, C in (2.29), as a function 
of q. Using Parseval’s relation (e.g. [21]) we have: 


(2.35) C < c 0 (T) sup 

*=7+trj, 


(* a + |*| a ) 1 ' a + ni 

s+(s> + \k\>)V' + \k\R(3/\k\) ' 


Here 7 > 0 is fixed and co(T) depends on 7 but is independent of the 
boundary condition. Dividing through by \k\ we see that our problem is 
reduced to the estimation of: 


(2.36) 


sup \Q(z)\, 

*(z)>0 


Q(z) = 


( z 2 + 1)1/2 + 1 
z + (z* + iy/* + R( z y 


Noting the form of R(z ), and recalling the choice of branch for the roots, 
we see that the Q is bounded as z —> 00 independent of q . Therefore, by 
the maximum principle, we can restrict attention to the imaginary axis, 
z = ir). In this case the poles of R are located on the imaginary axis with 
|t?| < 1 and \R\ is strictly decreasing in hi for hi > 1- Therefore, for hi > 1 
we have: 

IQI < — 1 + vV ~ 1 — < — ^ . 

Ini + vV - 1 - |fi(*)| 1 - WO I 

Now 



(2.38) 


1-I*(i)l = l 


2 

(9 + 1) 


9 



1 


To bound this we note: 


(2.39) 


7T 



dw 

Vl-w 2 


For q = 2p, even, we have: 


(2.40) 

Generally we have: 


' >2 E / 

j = 1 Jw 3 - 


dw 

2/U+l) Vl “ W 2 


(2.41) 


r w j 

Jwj-2 


dw 2 1 

> 


-2/( ? +l) VT-W> 9 + * y/l - ty 2 ’ 

and for j = 1, q > 2; 

(2.42) f" 7 == > ^-V^ I 2 1 

J~i VT^w* y/q + 1 9 + 1 \/l - tuj 
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Hence, for some constant, c, independent of q: 


(2.43) 


1 


< Cy/q. 


For q = 2pH- 1, odd, 

l-r* 

(2.44) ^ = 2 E/ „ , 


p+i 


dtu 

>/l — tu 3 ’ 


and (2.43) is similarly established. For |tj| < 1 we distinguish two cases, 
|rj| < |tuj.| and |??| > |tui|. For the former we have: 

(2.45) | < 2 | < - 7 ==j < cy/q, 

V 1 “ w l 

while for the latter it is easily shown that |Q| is maximized at |rj| = 1, 
which we have estimated above. Our final estimate is: 


(2.46) C < C(T)y/q, 

with C independent of q . 

To derive error estimates we first estimate the error in the quadrature 
formula for fixed t. Although the integrand doesn't possess a bounded 
derivative, its derivative is integrable. Therefore, employing the Peano 
Kernel representation of the error: 


(2.47) - (TTIh £ v /r -^ co ’”'‘ l £ 


Ci + 2 Cjt 


Here, the Cj are constants independent of t and q . Then we have: 

||£»r.p(t)||L 1 (( 0 ,T)) = J Q 1^7^ ~ {q+ l)x ^\/ 1_ cosw i t \ dt 

(2.48) < i( Cl T + c 2 T 2 ). 

Substituting this into (2.29) and using (2.46) we finally obtain, for some 
c(T) independent of q: 

(2.49) f ||c( , > t )il/f 1 [n3 < ^ - l » t )ll/r 3 [i^- 1 ]*- 

Jo 9 Jo 

We have shown, then, that the trapezoid rule produces a one- half order 
approximation to the true solution in the sense that the error decays at least 
as the square root of the reciprocal of the number of terms in the boundary 
condition. 
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A more accurate formula is obtained through the application of the 
Gaussian quadrature rule associated with the weight >/l — tt> 2 : 


(2.50) 


\(t) _ 1 - 2 j* * 3* 

& > sm cos Wit* w* = cos , 

t 9 + lf-r 9 + 1 1 3 9+1 


(2.51) 




1 V" ( sip2 ;+r) 2 


* 2 + w? 


This approximation is discussed in [15]. Remarkably, it is shown to 
be equivalent to the stable Pade approximants introduced by Engquist and 
Majda [7,8]. The stability of this approximation is well-known. To estimate 
the stability constant, we repeat the analysis given above for the trapezoid 
approximation. As in that case, we must estimate (1 — |R(t)|) _1 . We find: 

( ■i ? s [ n 2 

-* +1 - 

We also have: 


(2.53) 


Hence we conclude: 


^1 — w\ 


< c(g + 1). 


(2.54) 


C < C(T)q, 


with C independent of q. 

Using the error formula for the quadrature rule we find [15]: 


(2.55) 


IIW«)lk ( o,T)<e ^|^ -. 


From this, (2.29) and (2.54) we prove: 

( 2 - 56 ) ^ ^ e ('( 2 g + l )!) 3 J 0 

Here we see exponential convergence in q for smooth w , as with spec- 
tral approximations to the solutions of differential equations, so that the 
method might reasonably be called infinite order. 

Clearly, these estimates are nonuniform in T, as must be the case for 
spatially local boundary conditions. (See (2.18).) In [14,15] we presented a 
number of methods for approximating the convolution kernel by decaying 
functions, including exponential interpolation, exponential least squares 
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and approximations by Laguerre functions. None of these seemed entirely 
satisfactory: the interpolation could not be easily extended to high order, 
while the other approximations converged very slowly (if at all) in £i(0, oo). 
(We note that the proposed conditions have not yet been tested.) 

In order to use quadrature to construct long-time approximations we 
must derive an integral representation of Ji(t)/t involving an integrand 
which decays in t. This may be accomplished by treating (2.4) as a complex 
integral and deforming the contour. For example, let: 

(2.57) z = w + t(l - w 2 )P(w) 1 P > 0, w€[—l, 1]. 

Then: 

J -f = M/c (i “ 2l)1,v S 

(2.58) = — J \/l “ w 2 D(w,t)dw 1 

(2.59) D(w y t) = (A(ty)cos wt + / 2 (u>) sintut)e"^ 1 “ 1 " 3 ^ ti '^, 

(2.60) fi{w) = ff(u>) + (2 wP(w) — (1 — ™ 2 )P'(u/))/t(i{;), 

(2.61) /a(ii;) = (2tuP(u/) - (1 - w 2 )P'(w))g(w) - h(w) t 

(2.62) $(tu) = G(tu) + (1 - w 2 )P 2 (w)/G(w ), 


(2.63) G{w) = y/l + y/TTfTT ti>) 2 P 2 (tu)^l + v/l + (1 - w) a P 3 {w), 


(2.64) /»(tu) = (1 - tu)P(u>).ff(u>) - (1 + u>)P(u>)/JI(u>), 


(2.65) 


J 1 + v/TT (1 + w) 3 P 2 (wi 

H(w) = - 

y 1 + 4- (1 - ti») 2 P 2 (ti>) 


Since the function D(w> t) is smooth on the interval of integration, it 
is still reasonable to use the Gaussian quadrature scheme for the weight 
y/l — w 2 . This yields: 


iM I 




(2.66) 



11 


(2.67) 


sin 2 $r(/iK)(* + (1 - ™*)P{v>j)) + h(wj)wj) 
(H(i-«f)PW) J +^ 

We have not yet tested or analyzed these schemes. Finite time error 
estimates should be obtainable using the error formula for the quadrature 
rule, but the hope is that time independent estimates will hold. 

Problem B: Find a function or class of functions, P(iu), such that the 
resulting scheme leads to well-posed problems and kernels which converge 
to IC(t) in £i[(0, oo)]. 

3. Generalizations. 

3.1. Anisotropic Problems. The approximations discussed above 
are also applicable to problems with anisotropic wave propagation. As a 
first example, consider the convective wave equation with subsonic convec- 
tion: 

(S.!) ( u = V J U , Xl >L. 

We normalize so that = 1, 0 < M < 1. To formulate the exact 

boundary condition, we seek solutions of the form: 

(3.2) u= Ae Xx \ 

leading to the quadratic equation, 



*(*) = 


2(9 + !)^ 


£ r j( 2 ). 


(3.3) 


^5 4- iM = A 2 — |ife| 2 . 


The relevant solution, that is the one with negative real part for 9?(s) 
sufficiently large, is given by: 


A = 

(3.4) — (1 — Af 2 Wi) -1 ((1 -Mu l )s+ ((S 2 +(l-«V)|lf) 1/J -i)), 


5 = 5+ iM y ^ Lc>ikl. 


l> 1 


where 

(3.5) 
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To conveniently express approximations to this condition, we define the 
tangential material derivative^ ky 

(3.6) = -rr + ' ^tan» ^tan = A/ ( ( j >2 ) > 

La ut 

and let D/Dt denote the standard material derivative with respect to the 
full velocity U = M( u>i ... u> n ) T . The exact condition is then given by: 

(3.7) (1 + + i(l7tan • *)ti + (1 + Mui)\k\ 2 A * u, 

OX i ut 


_ Ji(\/1 - M 2 Ui\k\t 
y/l — 


Using the rational approximation y^l — ~Mhjf\k\R{z), with z — 5/{>J 1 — M 3 wJ|fc|), 
we obtain, in analogy with (2.11-2.12), 


(3.9) 


dv dv , 

(i + AfwO— + ^-+t(u t , 


P+1 

■*)* + E^* = 

3=1 


0, 


(3.10) + *(0un • k ) - V 1 - =Qj|fc| 2 (l + Mw 1 )v. 

The approximate condition is local in space under the same conditions for 
locality in the isotropic case, that is (2.13). Then we have: 


(3.11) 


dv Dv 

S>- =0, 

J = 1 


(3.12) (2k _ (1 - m 3 w 3 )/??v?„) = -(1 + v 3 „2^. 

We are confident that issues of consistency and convergence for these 
approximations could be handled as in the isotropic case, but we have not 
yet carried out the details. Below we see that the same operators appear 
in boundary conditions for the linearized compressible Euler and Navier- 
Stokes systems. 

3.2. Applications to Fluid Dynamics. We now consider the com- 
pressible Euler equations linearized about a uniform flow in two space di- 
mensions: 


Du dp 
~Dt + d x 


= 0, 


(3.13) 
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(3.14) 

Dv dp 
~Dt + dy ~ ’ 

(3.15) 

Dp du dv 
Di dx dy 


where D/Di = d/di -f Mw x d/dx + Mw v d/dy , u> 2 + a; 2 = 1, 0 < M < 1. 
We suppose u x > 0 and put artificial boundaries at x = ±L ; that is inflow 
at — L and outflow at L. 

To compute exact conditions we Laplace transform in t, Fourier trans- 
form in y, and rewrite the system in the form: 


< 3i < 


u 

v 

P 




Mw r s 

0 

s 


ik_ 


1-A f 7 ul 


MiJm 
ikMcj , 

1 -M*wl 


“l-jAuJ 

ik 

MWm 

Mw*.s 



Here, s = s + ikMu) y . Eigenvalues and left eigenvectors of the coefficient 
matrix are given by: 


j _ . _ A + Af w x s ,'t* * *. 

(3-17) Al= 1-M 2 ^ ’ 1=( * - A )« 

(3.18) A 2 = ~ j ~ > Ij=(* A), 

(3-19) A 3 = ^3 = 5 »*), 

where 

(3.20) i = (s 2 + (1 - M 2 * 2 )* 2 ) 1 / 2 . 

Noting that $ft(Ai) > 0 and < 0 for 3fc(s) > 0 we have, setting 

w — (u v p) T : 

(3.21) {[ w = 0, x = L\ = 1% w = 0, x — —L. 

The exact outflow boundary condition is then given by: 

(3.22) + ikMuty ^ (u — p) — ikMu) x v — (1 — M 2 u>l)k?A * p = 0. 

The exact inflow conditions are given by: 

(3.23) (— + (u-f p) — ikMw x v + (1 — M 2 w\)k 2 A * p = 0, 
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(3.24) 
Here, 

(3.25) 


(£ + iWfa, y ) 


v 4 ikMu} x u 4 ikp = 0. 


_ ikUUvt 
A{ ’ ) - M^lkt 


Using the rational approximation y'l — M 2 u}\\k\R{z) y for z = $/(%/! “ M 2 w* |fc|), 
equations (3.22- 3.23) become, in analogy with (3.9-3.10), 


( d \ 

— 4 ikMutyJ (u - p) - ilMwjti - = 0, 

(3.27) + tfcM Wy - Vl “ M^l\k\p^j hj = a^l - M 2 w 2 )p, 


at outflow and at inflow, 


(3.28) 


( 


— + ikMu 
dt 


y 



p+i 

— ikMu) x v 4 9j 
J'=i 


= 0 , 


(3.29) + ikMw y - y/T— M 2 w2|fc|pj) gj = a,-fc 2 ( 1 - M 2 « 2 )p. 

The approximate conditions are local in space under the same con- 
ditions for locality in the scalar case, that is (2.13). In particular, if we 
use the Engquist-Maj da-Pad e-Gaussian approximation, (2.51), we have, at 
outflow: 

( 3 - 30 ) _p) ~ Mui * - 5Z = °* 


(3.31) 




M 2 u)\ cos 




d 


q + \dy 


) = 


9 + 1 


* 2 /i »^2 2\^ 2 P 

sm 


The inflow conditions become: 


•Ptan 

~Dt 


dv 


(u + p) - Mw z ^ + Y^'l’j 

y i=i 


= 0, 


(3.32) 
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(3.33) 






M 2 l>1 cos 


-—7 sin 3 — — (1 - Af 3 w 3 ) 
+ 1 g + l v *' 


3* 
9 + 


t£) *< = 


,7& 

dy 2 ' 


(3.34) 


Dt 


du 

+ Mu x — + 


ap 


= 0. 


These conditions have been implemented by Goodrich [9] for channel 
flows, and shown to be very accurate even for moderate q. Note that 
(3.34) implies zero vorticity at inflow and, for the linearized problem, it is 
exact. Also, we have used an alternative representation of the approximate 
conditions with the property that the auxiliary functions are computed by 
solving first order equations. 

A similar construction has been carried out for the linearized, isen- 
tropic, compressible Navier- Stokes equations in the low Mach number limit 
by the author and Lorenz [18]. Here, we require two boundary conditions 
at outflow and three at inflow. The conditions (3.22), (3.23) and (3.34) 
remain the same to leading order. The additional condition at outflow is 


(3.35) 

and at inflow is given by: 


Dv dp 
~Dt + dy 


= 0 , 


(3.36) 


du dp 

&t + di 


= 0. 


We note that here we generally expect that long time computations, as mea- 
sured on the time scale of the sound waves, will be of interest. Therefore, 
nonlocal approximations may be efficient. As part of his doctoral disserta- 
tion, L. Xu is looking into such approximations and their applications in 
acoustics. Some of his results will appear in [19]. 

3.3. Corner Conditions. The boundary conditions, as discussed so 
far, only apply to half-space or periodic problems. To generalize their 
applicability, one must understand how to treat the case of an artificial 
boundary intersecting another part of the boundary at a corner. Collino 
[4] has solved this problem for the important special case of the isotropic 
wave equation, two artificial boundaries intersecting at a right angle, and 
spatially local boundary conditions. (We note that our formulation of the 
exact conditions and spatially nonlocal approximations is not valid in this 
case.) 

As a first attempt to generalize Collino’s results, we have considered 
the convective wave equation in a rectangular domain with all boundaries 
artificial and the spatially local boundary conditions (3.11-3.12). Collino’s 
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construction does not directly apply, as it is based on special exact solu- 
tions of the wave equation in a corner. We have, instead, looked at power 
series expansions in the corner. These do lead to compatibility conditions 
which c an be used to relate the auxiliary functions associated with distinct 
boundaries. However, it seems that the expansion must be carried out to 
order greater than 4 q to produce the q required conditions. We have tried 
to do this symbolically, but so far have succeeded only for q < 2. 

Problem C: Derive corner compatibility relations for local approximate 
boundary conditions for anisotropic systems. 

We also note that the theory of exact conditions and nonlocal approx- 
imations is still undeveloped. 

Problem D: Characterize the exact boundary condition for the wave equa- 
tion and convective wave equation with a rectangle or rectangular par- 
allelipiped as artificial boundary. Construct convergent temporally local 
approximate conditions. 

Another case of practical importance should be analyzed is the inter- 
section of an artificial boundary with a physical boundary, such as a solid 
wall. As a first example consider the convective wave equation in two 
space dimensions with convection in the x direction and walls at y = 0 and 
y = H: 

(3-37) (§t +M Jx) u=V2u ’ (*.y) € (O.oo) x (0,£T), 


(3.38) Q0-5 /3 o u = 0, y = 0; (»i— +/3ju=0, y=H. 

dy ay 

We can expand the solution: 

(3.39) u = 53 Yi (y) = Ai cos k * y + Bi sin 

i 

where the eigenvalues, Aj, are determined by the boundary conditions. 
Then, uj satisfies (3.7), and for the approximate local boundary condition 
(3.11)-(3.12) we have: 

(3.40) <t>j = hMy)- 

i 


Since the eigenfunctions, Yi, satisfy (3.38) so would the <f>j t under the as- 
sumption that the expansion is sufficiently regular. That is: 


d6i 

(3.4i) ~ = o- y = °; 


Qi-“ - <f>j = 0, 

dy 


y = S. 
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We should emphasize that this is generally not a true solution of the com- 
patibility problem, but simply an approximation which we have used. 

This approximation can be extended to the linearized, compressible 
Euler system, again with the physical assumption that the base flow is 
parallel to the walls, y = 0, S'. The physical boundary condition at this 
characteristic boundary is: 


(3.42) = 0, 

while a second relation, implied by the y-momentum equation is: 

dp 


(3.43) 


9y 


= 0 . 


Using these we see that the solution can be expanded in the form: 

u = ^ u\ cos Iny/H, v = Ez s ^ n 
i 

(3.44) p = Et pi cos foyfE- 

To relate these expansions to expansions of the auxiliary variables, we note 
that: 


(3.45) 


jt (q + 1 — i)x 

cos b cos — = 0. 

9 + 1 9 + 1 

Hence, equations (3.31) and (3.33) imply: 

(3.46) 4>j + <l> 9 +i-j = cos Jxy/jET, 


(3.47) Vi + V’j+i-i = 9 1 cos Ivy/ H, 

l 


with boundary conditions: 


(3.48) 


dy dy 1 dy dy 


These are the corner conditions used by Goodrich in [9]. 

3.4. Smooth Boundaries. Given the theoretical and practical dif- 
ficulties of high-order boundary conditions at artificial boundaries with 
corners, smooth artificial boundaries are an attractive alternative. Repre- 
sentations of the exact boundary condition for the wave equation at circular 
and spherical artificial boundaries are easily obtained. (See [17].) 

Using polar coordinates in two dimensions, we place the artificial bound- 
ary at r — R and expand the solution, u, in a Fourier series in 8: 


u(r,6,t) = ^ui e xli . 
\ 


(3.49) 
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The exact boundary condition is then given by: 

< 3 - 5l> ) lr + ^ + ^ a,+ 55 A(,/J1) * ii,=0 ' 

(3.51) A W = -*(»W +1 + s)- 

(Here, K\(z) is the modified Bessel function [!]•) 

We have not, as yet, found a closed form expression for A\ . We note 
that conditions based on far-field expansions (e.g. [2,16]) correspond to 
large z expansions of A\. By the results of [17], long time accuracy is 
difficult to achieve using time-local conditions, as Ao decays slowly as t — ► 
oc. 

In three dimensions, on the other hand, the exact condition takes a 
much simpler form. Here we use spherical coordinates and take r = R 
as our artificial boundary. The solution, u, is now expanded in spherical 
harmonics: 

(3.52) u(r,M,t) = 5Dii,y,(0,^), 

l 

where 

(3.53) ^«Ji = -l(l+l)Yi. 

In analogy with the two-dimensional case, the exact condition is related 
to the inverse Laplace transform of the logarithmic derivative of spherical 
Bessel functions. In particular we have: 

< 3 ' 54 > lr + ^ + 5 0, + 5i s,( ‘ / * ) * il, = 0 ' 


(3.55) 5i(«) 


( (z-^K l+1/2 (z)y i\ Pl (z) 

V z~ 1/2 Ki+ 1/aW *) Qi(*Y 


where, for / ^ 0, 


(3.56) Pt{z) = 


t-i 


E 


(21 - Jfc)! 
Jk!(J — Jb — 1)! 


(2*)\ 


i 


<?«(*) = E 

k= o 


(2/ - Jb)! 
Jb!(l — *)! 


(2 z ) k i 


and So = 0. 

We immediately observe that Sj(z ) is rational and, therefore, (3.54) 
may be localized in time. (The fact that the exact boundary condition 
can be localized in time for a finite number of spherical harmonics is also 
noted by Grote and Keller [10].) What has not been accomplished so 
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far is to derive a convenient factorization or expansion of Pi{z)/Qi{z) } to 
facilitate numerical implementation. However, even if this must be done 
numerically, the availability of simple, exact conditions makes the use of a 
spherical boundary very attractive. 

We note that an important advantage of rectangular boundaries is 
the possibility of adjusting the aspect ratio of the domain. To do this 
with a smooth boundary, while maintaining a natural, separable coordinate 
system, one might need elliptical and spheroidal coordinates. This suggests 
our final problem. 

Problem E: Characterize and approximate the exact boundary condition 
for the wave equation at elliptical and spheroidal boundaries. 
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A. Derivation of Equation (2.4). Our goal is to compute a useful 
representation of the inverse Laplace transform of: 

(A.1) £(,)= (,2 + i) 1/2 -,. 

Although the inverse is available in standard tables, its importance to our 
work suggests the inclusion of a direct verification. 

Following the ideas given in [5, Ch. 38], we will derive a simple differ- 
ential equation satisfied by /C(t). We first note that K satisfies a first order 
differential equation with quadratic coefficients: 

(A.2) (’’ + l)f =»*-!• 

Recall that: 

(A.3) ^ = W). 

Therefore, (A.2) is equivalent to: 

(A. 4) ^(tK)+^+tK = l-2fC(0 + ). 

That is, fC(t) satisfies the differential equation: 

(A.5) t^ + 3^-+tfC = MK = 0, 

subject to the initial condition, 

( a . 6 ) m=\> 

in addition to a growth condition at infinity. 



We now verify that: 


1 

(A.7) K(t) = - / VT — w 2 cos wtdw, 

* J - 1 

solves this problem. That the initial condition is satisfied may be checked 
directly. Applying the differential operator we obtain: 


■kMK 

= f (1 — w 2 ) 3 / 2 t cos wtdw — 3 

1 yjl — w 2 w sin wtdw 


J- 1 

J-i 

(A.8) 

— J (1 — w 2 ) 3 t 2 t cos wtdw + j 

-p-(l — w 2 ) 3 / 2 sin wtdw 
dw 


= 0 . 


Finally, we 

recall the identity, 


(A.9) 

Ji(t) ~ ~ J \/l — w 2 ( 

:os wtdw , 


which may be found, e.g., in [1]. 
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